Set working directory

Visualizations

library(tidyverse)
## ── Attaching packages ────────────── tidyverse 1.2.1 ──
## ✔ ggplot2 3.0.0     ✔ purrr   0.2.5
## ✔ tibble  1.4.2     ✔ dplyr   0.7.6
## ✔ tidyr   0.8.1     ✔ stringr 1.3.1
## ✔ readr   1.1.1     ✔ forcats 0.3.0
## ── Conflicts ───────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
library(ggthemes)

plant_data <- read.csv('plant_data.csv', header = TRUE)

ggplot(plant_data, 
       aes(x=Biostimulant,
           y=Leaves/1000))+ #,fill=TimePoint)) +
  facet_grid(.~TimePoint) +
  geom_violin() +
  xlab('Biostimulants') +
  ylab('Leaves dry weight (g)') + theme_igray() +
  theme(axis.text.x  = element_text(angle=90, vjust=0.5)) +
  scale_fill_pander()
## Warning: Removed 1 rows containing non-finite values (stat_ydensity).

Data Cleaning

##Check NANs
colSums(is.na(plant_data))
##   X.SampleID    TimePoint Biostimulant          Rep         Area 
##            0            0            0            0            1 
##     SPAD.avg       Leaves        Stems         Root          Tot 
##            0            1            0            0            0
plant_data[!complete.cases(plant_data),]

Imputing missing data

##Imputing NANs with averages of subsets of the data (based on timepoint and biostimulants, thus taking the mean from replicates of the same missing sample)

plant_data$Area <- ifelse(is.na(plant_data$Area),
                      ave(plant_data[(plant_data$TimePoint==plant_data[which(is.na(plant_data$Area)),'TimePoint'] & plant_data$Biostimulant==plant_data[which(is.na(plant_data$Area)),'Biostimulant']),'Area'],
                                 FUN = function(x) mean(x, na.rm = TRUE)),plant_data$Area)

plant_data$Leaves <- ifelse(is.na(plant_data$Leaves),
                      ave(plant_data[(plant_data$TimePoint==plant_data[which(is.na(plant_data$Leaves)),'TimePoint'] & plant_data$Biostimulant==plant_data[which(is.na(plant_data$Leaves)),'Biostimulant']),'Leaves'],
                                 FUN = function(x) mean(x, na.rm = TRUE)),plant_data$Leaves)

##Check NANs
colSums(is.na(plant_data))
##   X.SampleID    TimePoint Biostimulant          Rep         Area 
##            0            0            0            0            0 
##     SPAD.avg       Leaves        Stems         Root          Tot 
##            0            0            0            0            0

Feature engineering

Combine categorical variables

lda_variables <- plant_data %>%
mutate(Biostimulant_TimePoint=paste(Biostimulant, TimePoint, sep = '_'))

LDA Area vs TimePoint

library(MASS)
## 
## Attaching package: 'MASS'
## The following object is masked from 'package:dplyr':
## 
##     select
area_vs_treatment_set = dplyr::select(plant_data, Area, TimePoint)
# Feature Scaling
area_vs_treatment_set[1] = scale(area_vs_treatment_set[1])
# Applying LDA
lda_area_vs_treatment = lda(formula = TimePoint ~ Area, 
                              data = area_vs_treatment_set)
area_vs_treatment_results = as.data.frame(predict(lda_area_vs_treatment, area_vs_treatment_set))
area_vs_treatment_results = area_vs_treatment_results[c(6, 1)]
names(area_vs_treatment_results) <- c('LDA1_area_vs_treatment')

# Check how much variance the LDA covers
lda_area_vs_treatment$svd^2/sum(lda_area_vs_treatment$svd^2)
## [1] 1
lda_variables <- cbind(lda_variables, area_vs_treatment_results[1])

LDA Biomass vs Timepoint

biomass_vs_treatment_set = dplyr::select(plant_data, Leaves, Root, Stems, TimePoint)
# Feature Scaling
biomass_vs_treatment_set[-4] = scale(biomass_vs_treatment_set[-4])
# Applying LDA
lda_biomass_vs_treatment = lda(formula = TimePoint ~ Leaves+Root+Stems, 
                              data = biomass_vs_treatment_set)
biomass_vs_treatment_results = as.data.frame(predict(lda_biomass_vs_treatment, biomass_vs_treatment_set))
biomass_vs_treatment_results = biomass_vs_treatment_results[c(6,7,8,1)]
names(biomass_vs_treatment_results) <- c('LDA1_biomass_vs_treatment')
lda_biomass_vs_treatment$svd^2/sum(lda_biomass_vs_treatment$svd^2)
## [1] 0.972823981 0.025374543 0.001801476
lda_variables <- cbind(lda_variables, biomass_vs_treatment_results[1])

Statistical exploratory analysis

Correlation plots

library(corrplot)
## corrplot 0.84 loaded
dummyvars <- model.matrix(X.SampleID~Biostimulant+Biostimulant_TimePoint,lda_variables)
encoded_vars <- cbind(dummyvars[,-1], lda_variables[,which(names(lda_variables)!='LinkerPrimerSequence' &
                                                  names(lda_variables)!='BarcodeSequence' &
                                                  names(lda_variables)!='Description' &
                                                  names(lda_variables)!='InputFileName' &
                                                  names(lda_variables)!='X.SampleID' &
                                                  names(lda_variables)!='Rep' &
                                                  names(lda_variables)!='Biostimulant' &
                                                  names(lda_variables)!='TimePoint' &  
                                                  names(lda_variables)!='Biostimulant_TimePoint')])
corr.mat = cor(encoded_vars, use='complete.obs')
p.mat = cor.mtest(encoded_vars)$p
row.names(p.mat) <- row.names(corr.mat)
colnames(p.mat) <- colnames(corr.mat)
corrplot(corr.mat, type='lower', 
         tl.pos = 'lt', 
         addCoefasPercent = TRUE, 
         method = 'square',
         number.cex = .4, tl.cex = 0.8,
         order = 'FPC',  tl.col="black",
         p.mat = p.mat, sig.level = 0.01, 
         insig = 'label_sig', 
         pch.cex = .3, 
         pch = '*') +
corrplot(p.mat, type = 'upper', add = T, tl.pos = 'n', method = 'circle', 
         col =rainbow(40, start = 5/6, end = 4/6), cl.lim =c(0,1), diag = FALSE)

##                                        Leaves          Tot        Root
## Leaves                             1.00000000  0.995723837  0.96940097
## Tot                                0.99572384  1.000000000  0.99375599
## Root                               0.96940097  0.993755990  1.00000000
## Stems                              0.98362181  1.002198386  0.98875469
## LDA1_biomass_vs_treatment          1.16468160  1.193950825  1.21145877
## Area                               0.93359382  0.864161010  0.81455220
## LDA1_area_vs_treatment             0.93359382  0.864161010  0.81455220
## Biostimulant_TimePointSumagrow_T3  0.34340703  0.586056253  0.60882466
## Biostimulant_TimePointControl_T3   0.39367886  0.611397006  0.59845048
## Biostimulant_TimePointInocucor_T3  0.35351495  0.580923657  0.58575973
## Biostimulant_TimePointPathway_T3   0.34334009  0.585020527  0.58565901
## Biostimulant_TimePointEndomaxx_T3  0.51651659  0.324622323  0.58424893
## SPAD.avg                           0.43433526  0.155594941  0.35266336
## Biostimulant_TimePointSumagrow_T2  0.24553066 -0.003210153  0.19821858
## Biostimulant_TimePointPathway_T2   0.24093189  0.006330031  0.21500864
## BiostimulantSumagrow               0.22194174  0.225156785  0.02630018
## Biostimulant_TimePointInocucor_T2  0.25237861  0.200961960 -0.06869508
## Biostimulant_TimePointControl_T2   0.22641263  0.178532746 -0.08003617
## BiostimulantPathway                0.21122233  0.225802627  0.02099342
## BiostimulantInocucor               0.23159250  0.218743512  0.19786498
## Biostimulant_TimePointEndomaxx_T2  0.19426812  0.182758507  0.16405429
## BiostimulantEndomaxx               0.16011999  0.175250132  0.20207966
## Biostimulant_TimePointInocucor_T1  0.09170286  0.078345856  0.08324646
## Biostimulant_TimePointSumagrow_T1  0.69918270  0.625862955  0.67662286
## Biostimulant_TimePointPathway_T1  -0.06037137  0.450450390  0.30431481
## Biostimulant_TimePointControl_T1   0.44454594  0.686100993  0.84565505
## Biostimulant_TimePointEndomaxx_T1  0.53233412  0.787571033  0.68139120
## Biostimulant_TimePointPathway_T0   0.67526442  0.665134628  0.61804883
## Biostimulant_TimePointSumagrow_T0  0.45171076  0.732569343  0.68105724
## Biostimulant_TimePointInocucor_T0  0.58350393  0.550206651  0.60703526
## Biostimulant_TimePointEndomaxx_T0  0.72625469  0.624205873  0.59059149
##                                          Stems LDA1_biomass_vs_treatment
## Leaves                             0.983621806                1.16468160
## Tot                                1.002198386                1.19395082
## Root                               0.988754694                1.21145877
## Stems                              1.000000000                1.19011866
## LDA1_biomass_vs_treatment          1.190118657                1.00000000
## Area                               0.823698571                1.14925044
## LDA1_area_vs_treatment             0.823698571                1.14925044
## Biostimulant_TimePointSumagrow_T3  0.591388340                0.97002406
## Biostimulant_TimePointControl_T3   0.615330622                0.95106366
## Biostimulant_TimePointInocucor_T3  0.583724695                0.94233301
## Biostimulant_TimePointPathway_T3   0.601218946                0.94390155
## Biostimulant_TimePointEndomaxx_T3  0.522098423                0.94730194
## SPAD.avg                           0.320624190                0.69525414
## Biostimulant_TimePointSumagrow_T2  0.185812419                0.54702235
## Biostimulant_TimePointPathway_T2   0.202540305                0.56688202
## BiostimulantSumagrow               0.220367806                0.59647682
## Biostimulant_TimePointInocucor_T2  0.190674869                0.48410639
## Biostimulant_TimePointControl_T2   0.165622680                0.47486454
## BiostimulantPathway                0.233099163                0.59218386
## BiostimulantInocucor               0.007143058                0.54902220
## Biostimulant_TimePointEndomaxx_T2 -0.028229295                0.51591908
## BiostimulantEndomaxx              -0.038512433                0.56393470
## Biostimulant_TimePointInocucor_T1 -0.144610985                0.43846561
## Biostimulant_TimePointSumagrow_T1  0.556747101                0.77190212
## Biostimulant_TimePointPathway_T1  -0.022550517               -0.10627022
## Biostimulant_TimePointControl_T1   0.755480590                0.01715083
## Biostimulant_TimePointEndomaxx_T1  0.770097755               -0.05805714
## Biostimulant_TimePointPathway_T0   0.583524966               -0.08345161
## Biostimulant_TimePointSumagrow_T0  0.696905277               -0.08708204
## Biostimulant_TimePointInocucor_T0  0.510865654                0.70979950
## Biostimulant_TimePointEndomaxx_T0  0.565031908               -0.08952663
##                                          Area LDA1_area_vs_treatment
## Leaves                             0.93359382             0.93359382
## Tot                                0.86416101             0.86416101
## Root                               0.81455220             0.81455220
## Stems                              0.82369857             0.82369857
## LDA1_biomass_vs_treatment          1.14925044             1.14925044
## Area                               1.00000000             1.56805979
## LDA1_area_vs_treatment             1.56805979             1.00000000
## Biostimulant_TimePointSumagrow_T3  0.78449689             0.78449689
## Biostimulant_TimePointControl_T3   0.68032365             0.68032365
## Biostimulant_TimePointInocucor_T3  0.71707683             0.71707683
## Biostimulant_TimePointPathway_T3   0.71058866             0.71058866
## Biostimulant_TimePointEndomaxx_T3  0.76352779             0.76352779
## SPAD.avg                           1.11187992             1.11187992
## Biostimulant_TimePointSumagrow_T2  0.72931728             0.72931728
## Biostimulant_TimePointPathway_T2   0.71262026             0.71262026
## BiostimulantSumagrow               0.60320515             0.60320515
## Biostimulant_TimePointInocucor_T2  0.72457336             0.72457336
## Biostimulant_TimePointControl_T2   0.73897053             0.73897053
## BiostimulantPathway                0.54579396             0.54579396
## BiostimulantInocucor               0.59601466             0.59601466
## Biostimulant_TimePointEndomaxx_T2  0.69574722             0.69574722
## BiostimulantEndomaxx               0.54720190             0.54720190
## Biostimulant_TimePointInocucor_T1  0.63552361             0.63552361
## Biostimulant_TimePointSumagrow_T1  0.08381664             0.24399966
## Biostimulant_TimePointPathway_T1   0.13260546             0.30249675
## Biostimulant_TimePointControl_T1   0.88853235             0.01096427
## Biostimulant_TimePointEndomaxx_T1  0.58601563            -0.02740427
## Biostimulant_TimePointPathway_T0   0.04483800            -0.34001125
## Biostimulant_TimePointSumagrow_T0  0.38003179            -0.33501550
## Biostimulant_TimePointInocucor_T0 -0.25970187            -0.09951884
## Biostimulant_TimePointEndomaxx_T0 -0.02269113            -0.33401142
##                                   Biostimulant_TimePointSumagrow_T3
## Leaves                                                   0.34340703
## Tot                                                      0.58605625
## Root                                                     0.60882466
## Stems                                                    0.59138834
## LDA1_biomass_vs_treatment                                0.97002406
## Area                                                     0.78449689
## LDA1_area_vs_treatment                                   0.78449689
## Biostimulant_TimePointSumagrow_T3                        1.00000000
## Biostimulant_TimePointControl_T3                         0.51542821
## Biostimulant_TimePointInocucor_T3                        0.51542821
## Biostimulant_TimePointPathway_T3                         0.51542821
## Biostimulant_TimePointEndomaxx_T3                        0.51542821
## SPAD.avg                                                 0.71247714
## Biostimulant_TimePointSumagrow_T2                        0.51542821
## Biostimulant_TimePointPathway_T2                         0.51542821
## BiostimulantSumagrow                                     1.02689126
## Biostimulant_TimePointInocucor_T2                        0.51542821
## Biostimulant_TimePointControl_T2                         0.51542821
## BiostimulantPathway                                      0.45335192
## BiostimulantInocucor                                     0.45335192
## Biostimulant_TimePointEndomaxx_T2                        0.51542821
## BiostimulantEndomaxx                                     0.45335192
## Biostimulant_TimePointInocucor_T1                        0.51542821
## Biostimulant_TimePointSumagrow_T1                       -0.05243893
## Biostimulant_TimePointPathway_T1                        -0.05261428
## Biostimulant_TimePointControl_T1                        -0.04379979
## Biostimulant_TimePointEndomaxx_T1                       -0.01792098
## Biostimulant_TimePointPathway_T0                        -0.02652675
## Biostimulant_TimePointSumagrow_T0                       -0.03270557
## Biostimulant_TimePointInocucor_T0                       -0.05243893
## Biostimulant_TimePointEndomaxx_T0                       -0.02123148
##                                   Biostimulant_TimePointControl_T3
## Leaves                                                  0.39367886
## Tot                                                     0.61139701
## Root                                                    0.59845048
## Stems                                                   0.61533062
## LDA1_biomass_vs_treatment                               0.95106366
## Area                                                    0.68032365
## LDA1_area_vs_treatment                                  0.68032365
## Biostimulant_TimePointSumagrow_T3                       0.51542821
## Biostimulant_TimePointControl_T3                        1.00000000
## Biostimulant_TimePointInocucor_T3                       0.51542821
## Biostimulant_TimePointPathway_T3                        0.51542821
## Biostimulant_TimePointEndomaxx_T3                       0.51542821
## SPAD.avg                                                0.47076315
## Biostimulant_TimePointSumagrow_T2                       0.51542821
## Biostimulant_TimePointPathway_T2                        0.51542821
## BiostimulantSumagrow                                    0.45335192
## Biostimulant_TimePointInocucor_T2                       0.51542821
## Biostimulant_TimePointControl_T2                        0.51542821
## BiostimulantPathway                                     0.45335192
## BiostimulantInocucor                                    0.45335192
## Biostimulant_TimePointEndomaxx_T2                       0.51542821
## BiostimulantEndomaxx                                    0.45335192
## Biostimulant_TimePointInocucor_T1                       0.51542821
## Biostimulant_TimePointSumagrow_T1                       0.71370780
## Biostimulant_TimePointPathway_T1                        0.01843060
## Biostimulant_TimePointControl_T1                        0.06427478
## Biostimulant_TimePointEndomaxx_T1                       0.02945594
## Biostimulant_TimePointPathway_T0                        0.07628471
## Biostimulant_TimePointSumagrow_T0                       0.04773485
## Biostimulant_TimePointInocucor_T0                       0.71370780
## Biostimulant_TimePointEndomaxx_T0                       0.07993096
##                                   Biostimulant_TimePointInocucor_T3
## Leaves                                                    0.3535149
## Tot                                                       0.5809237
## Root                                                      0.5857597
## Stems                                                     0.5837247
## LDA1_biomass_vs_treatment                                 0.9423330
## Area                                                      0.7170768
## LDA1_area_vs_treatment                                    0.7170768
## Biostimulant_TimePointSumagrow_T3                         0.5154282
## Biostimulant_TimePointControl_T3                          0.5154282
## Biostimulant_TimePointInocucor_T3                         1.0000000
## Biostimulant_TimePointPathway_T3                          0.5154282
## Biostimulant_TimePointEndomaxx_T3                         0.5154282
## SPAD.avg                                                  0.5783072
## Biostimulant_TimePointSumagrow_T2                         0.5154282
## Biostimulant_TimePointPathway_T2                          0.5154282
## BiostimulantSumagrow                                      0.4533519
## Biostimulant_TimePointInocucor_T2                         0.5154282
## Biostimulant_TimePointControl_T2                          0.5154282
## BiostimulantPathway                                       0.4533519
## BiostimulantInocucor                                      1.0268913
## Biostimulant_TimePointEndomaxx_T2                         0.5154282
## BiostimulantEndomaxx                                      0.4533519
## Biostimulant_TimePointInocucor_T1                         0.5154282
## Biostimulant_TimePointSumagrow_T1                         0.1119577
## Biostimulant_TimePointPathway_T1                          0.5086011
## Biostimulant_TimePointControl_T1                          0.7932143
## Biostimulant_TimePointEndomaxx_T1                         0.7069257
## Biostimulant_TimePointPathway_T0                          0.5488853
## Biostimulant_TimePointSumagrow_T0                         0.6969154
## Biostimulant_TimePointInocucor_T0                         0.1119577
## Biostimulant_TimePointEndomaxx_T0                         0.5190482
##                                   Biostimulant_TimePointPathway_T3
## Leaves                                                  0.34334009
## Tot                                                     0.58502053
## Root                                                    0.58565901
## Stems                                                   0.60121895
## LDA1_biomass_vs_treatment                               0.94390155
## Area                                                    0.71058866
## LDA1_area_vs_treatment                                  0.71058866
## Biostimulant_TimePointSumagrow_T3                       0.51542821
## Biostimulant_TimePointControl_T3                        0.51542821
## Biostimulant_TimePointInocucor_T3                       0.51542821
## Biostimulant_TimePointPathway_T3                        1.00000000
## Biostimulant_TimePointEndomaxx_T3                       0.51542821
## SPAD.avg                                                0.53457190
## Biostimulant_TimePointSumagrow_T2                       0.51542821
## Biostimulant_TimePointPathway_T2                        0.51542821
## BiostimulantSumagrow                                    0.45335192
## Biostimulant_TimePointInocucor_T2                       0.51542821
## Biostimulant_TimePointControl_T2                        0.51542821
## BiostimulantPathway                                     1.02689126
## BiostimulantInocucor                                    0.45335192
## Biostimulant_TimePointEndomaxx_T2                       0.51542821
## BiostimulantEndomaxx                                    0.45335192
## Biostimulant_TimePointInocucor_T1                       0.51542821
## Biostimulant_TimePointSumagrow_T1                      -0.02023767
## Biostimulant_TimePointPathway_T1                        0.70258484
## Biostimulant_TimePointControl_T1                       -0.05190465
## Biostimulant_TimePointEndomaxx_T1                      -0.05205926
## Biostimulant_TimePointPathway_T0                       -0.05260295
## Biostimulant_TimePointSumagrow_T0                      -0.05233390
## Biostimulant_TimePointInocucor_T0                      -0.02023767
## Biostimulant_TimePointEndomaxx_T0                      -0.05261217
##                                   Biostimulant_TimePointEndomaxx_T3
## Leaves                                                   0.51651659
## Tot                                                      0.32462232
## Root                                                     0.58424893
## Stems                                                    0.52209842
## LDA1_biomass_vs_treatment                                0.94730194
## Area                                                     0.76352779
## LDA1_area_vs_treatment                                   0.76352779
## Biostimulant_TimePointSumagrow_T3                        0.51542821
## Biostimulant_TimePointControl_T3                         0.51542821
## Biostimulant_TimePointInocucor_T3                        0.51542821
## Biostimulant_TimePointPathway_T3                         0.51542821
## Biostimulant_TimePointEndomaxx_T3                        1.00000000
## SPAD.avg                                                 0.53930425
## Biostimulant_TimePointSumagrow_T2                        0.51542821
## Biostimulant_TimePointPathway_T2                         0.51542821
## BiostimulantSumagrow                                     0.45335192
## Biostimulant_TimePointInocucor_T2                        0.51542821
## Biostimulant_TimePointControl_T2                         0.51542821
## BiostimulantPathway                                      0.45335192
## BiostimulantInocucor                                     0.45335192
## Biostimulant_TimePointEndomaxx_T2                        0.51542821
## BiostimulantEndomaxx                                     1.02689126
## Biostimulant_TimePointInocucor_T1                        0.51542821
## Biostimulant_TimePointSumagrow_T1                       -0.05229159
## Biostimulant_TimePointPathway_T1                        -0.05220574
## Biostimulant_TimePointControl_T1                        -0.04365947
## Biostimulant_TimePointEndomaxx_T1                       -0.01724260
## Biostimulant_TimePointPathway_T0                        -0.02586869
## Biostimulant_TimePointSumagrow_T0                       -0.03230415
## Biostimulant_TimePointInocucor_T0                       -0.05229159
## Biostimulant_TimePointEndomaxx_T0                       -0.02041589
##                                       SPAD.avg
## Leaves                             0.434335261
## Tot                                0.155594941
## Root                               0.352663355
## Stems                              0.320624190
## LDA1_biomass_vs_treatment          0.695254139
## Area                               1.111879922
## LDA1_area_vs_treatment             1.111879922
## Biostimulant_TimePointSumagrow_T3  0.712477142
## Biostimulant_TimePointControl_T3   0.470763146
## Biostimulant_TimePointInocucor_T3  0.578307153
## Biostimulant_TimePointPathway_T3   0.534571897
## Biostimulant_TimePointEndomaxx_T3  0.539304249
## SPAD.avg                           1.000000000
## Biostimulant_TimePointSumagrow_T2  0.776597916
## Biostimulant_TimePointPathway_T2   0.691311566
## BiostimulantSumagrow               0.708351808
## Biostimulant_TimePointInocucor_T2  0.597704597
## Biostimulant_TimePointControl_T2   0.710552999
## BiostimulantPathway                0.639186188
## BiostimulantInocucor               0.518819309
## Biostimulant_TimePointEndomaxx_T2  0.514498402
## BiostimulantEndomaxx               0.405536168
## Biostimulant_TimePointInocucor_T1  0.754444267
## Biostimulant_TimePointSumagrow_T1  0.638669744
## Biostimulant_TimePointPathway_T1   0.233784243
## Biostimulant_TimePointControl_T1   0.389563933
## Biostimulant_TimePointEndomaxx_T1  0.280435362
## Biostimulant_TimePointPathway_T0   0.008915816
## Biostimulant_TimePointSumagrow_T0 -0.125088022
## Biostimulant_TimePointInocucor_T0  0.147441171
## Biostimulant_TimePointEndomaxx_T0 -0.223017694
##                                   Biostimulant_TimePointSumagrow_T2
## Leaves                                                  0.245530658
## Tot                                                    -0.003210153
## Root                                                    0.198218581
## Stems                                                   0.185812419
## LDA1_biomass_vs_treatment                               0.547022353
## Area                                                    0.729317284
## LDA1_area_vs_treatment                                  0.729317284
## Biostimulant_TimePointSumagrow_T3                       0.515428211
## Biostimulant_TimePointControl_T3                        0.515428211
## Biostimulant_TimePointInocucor_T3                       0.515428211
## Biostimulant_TimePointPathway_T3                        0.515428211
## Biostimulant_TimePointEndomaxx_T3                       0.515428211
## SPAD.avg                                                0.776597916
## Biostimulant_TimePointSumagrow_T2                       1.000000000
## Biostimulant_TimePointPathway_T2                        0.515428211
## BiostimulantSumagrow                                    1.026891258
## Biostimulant_TimePointInocucor_T2                       0.515428211
## Biostimulant_TimePointControl_T2                        0.515428211
## BiostimulantPathway                                     0.453351923
## BiostimulantInocucor                                    0.453351923
## Biostimulant_TimePointEndomaxx_T2                       0.515428211
## BiostimulantEndomaxx                                    0.453351923
## Biostimulant_TimePointInocucor_T1                       0.515428211
## Biostimulant_TimePointSumagrow_T1                       0.035169636
## Biostimulant_TimePointPathway_T1                        0.695264952
## Biostimulant_TimePointControl_T1                        0.610440398
## Biostimulant_TimePointEndomaxx_T1                       0.762832520
## Biostimulant_TimePointPathway_T0                        0.403323976
## Biostimulant_TimePointSumagrow_T0                       0.850386420
## Biostimulant_TimePointInocucor_T0                       0.035169636
## Biostimulant_TimePointEndomaxx_T0                       0.309322388
##                                   Biostimulant_TimePointPathway_T2
## Leaves                                                 0.240931890
## Tot                                                    0.006330031
## Root                                                   0.215008640
## Stems                                                  0.202540305
## LDA1_biomass_vs_treatment                              0.566882017
## Area                                                   0.712620257
## LDA1_area_vs_treatment                                 0.712620257
## Biostimulant_TimePointSumagrow_T3                      0.515428211
## Biostimulant_TimePointControl_T3                       0.515428211
## Biostimulant_TimePointInocucor_T3                      0.515428211
## Biostimulant_TimePointPathway_T3                       0.515428211
## Biostimulant_TimePointEndomaxx_T3                      0.515428211
## SPAD.avg                                               0.691311566
## Biostimulant_TimePointSumagrow_T2                      0.515428211
## Biostimulant_TimePointPathway_T2                       1.000000000
## BiostimulantSumagrow                                   0.453351923
## Biostimulant_TimePointInocucor_T2                      0.515428211
## Biostimulant_TimePointControl_T2                       0.515428211
## BiostimulantPathway                                    1.026891258
## BiostimulantInocucor                                   0.453351923
## Biostimulant_TimePointEndomaxx_T2                      0.515428211
## BiostimulantEndomaxx                                   0.453351923
## Biostimulant_TimePointInocucor_T1                      0.515428211
## Biostimulant_TimePointSumagrow_T1                      0.051659546
## Biostimulant_TimePointPathway_T1                       0.858919594
## Biostimulant_TimePointControl_T1                      -0.052556849
## Biostimulant_TimePointEndomaxx_T1                     -0.052602140
## Biostimulant_TimePointPathway_T0                      -0.052605175
## Biostimulant_TimePointSumagrow_T0                     -0.052597425
## Biostimulant_TimePointInocucor_T0                      0.051659546
## Biostimulant_TimePointEndomaxx_T0                     -0.052606172
##                                   BiostimulantSumagrow
## Leaves                                      0.22194174
## Tot                                         0.22515678
## Root                                        0.02630018
## Stems                                       0.22036781
## LDA1_biomass_vs_treatment                   0.59647682
## Area                                        0.60320515
## LDA1_area_vs_treatment                      0.60320515
## Biostimulant_TimePointSumagrow_T3           1.02689126
## Biostimulant_TimePointControl_T3            0.45335192
## Biostimulant_TimePointInocucor_T3           0.45335192
## Biostimulant_TimePointPathway_T3            0.45335192
## Biostimulant_TimePointEndomaxx_T3           0.45335192
## SPAD.avg                                    0.70835181
## Biostimulant_TimePointSumagrow_T2           1.02689126
## Biostimulant_TimePointPathway_T2            0.45335192
## BiostimulantSumagrow                        1.00000000
## Biostimulant_TimePointInocucor_T2           0.45335192
## Biostimulant_TimePointControl_T2            0.45335192
## BiostimulantPathway                         0.31805979
## BiostimulantInocucor                        0.31805979
## Biostimulant_TimePointEndomaxx_T2           0.45335192
## BiostimulantEndomaxx                        0.31805979
## Biostimulant_TimePointInocucor_T1           0.45335192
## Biostimulant_TimePointSumagrow_T1           0.45897637
## Biostimulant_TimePointPathway_T1           -0.01609344
## Biostimulant_TimePointControl_T1           -0.10566186
## Biostimulant_TimePointEndomaxx_T1          -0.07906615
## Biostimulant_TimePointPathway_T0           -0.08771212
## Biostimulant_TimePointSumagrow_T0           0.47932140
## Biostimulant_TimePointInocucor_T0          -0.11456297
## Biostimulant_TimePointEndomaxx_T0          -0.08221434
##                                   Biostimulant_TimePointInocucor_T2
## Leaves                                                   0.25237861
## Tot                                                      0.20096196
## Root                                                    -0.06869508
## Stems                                                    0.19067487
## LDA1_biomass_vs_treatment                                0.48410639
## Area                                                     0.72457336
## LDA1_area_vs_treatment                                   0.72457336
## Biostimulant_TimePointSumagrow_T3                        0.51542821
## Biostimulant_TimePointControl_T3                         0.51542821
## Biostimulant_TimePointInocucor_T3                        0.51542821
## Biostimulant_TimePointPathway_T3                         0.51542821
## Biostimulant_TimePointEndomaxx_T3                        0.51542821
## SPAD.avg                                                 0.59770460
## Biostimulant_TimePointSumagrow_T2                        0.51542821
## Biostimulant_TimePointPathway_T2                         0.51542821
## BiostimulantSumagrow                                     0.45335192
## Biostimulant_TimePointInocucor_T2                        1.00000000
## Biostimulant_TimePointControl_T2                         0.51542821
## BiostimulantPathway                                      0.45335192
## BiostimulantInocucor                                     1.02689126
## Biostimulant_TimePointEndomaxx_T2                        0.51542821
## BiostimulantEndomaxx                                     0.45335192
## Biostimulant_TimePointInocucor_T1                        0.51542821
## Biostimulant_TimePointSumagrow_T1                        0.84326922
## Biostimulant_TimePointPathway_T1                        -0.01723015
## Biostimulant_TimePointControl_T1                         0.08445510
## Biostimulant_TimePointEndomaxx_T1                        0.05172942
## Biostimulant_TimePointPathway_T0                         0.08648990
## Biostimulant_TimePointSumagrow_T0                        0.06700611
## Biostimulant_TimePointInocucor_T0                        0.84326922
## Biostimulant_TimePointEndomaxx_T0                        0.08880906
##                                   Biostimulant_TimePointControl_T2
## Leaves                                                  0.22641263
## Tot                                                     0.17853275
## Root                                                   -0.08003617
## Stems                                                   0.16562268
## LDA1_biomass_vs_treatment                               0.47486454
## Area                                                    0.73897053
## LDA1_area_vs_treatment                                  0.73897053
## Biostimulant_TimePointSumagrow_T3                       0.51542821
## Biostimulant_TimePointControl_T3                        0.51542821
## Biostimulant_TimePointInocucor_T3                       0.51542821
## Biostimulant_TimePointPathway_T3                        0.51542821
## Biostimulant_TimePointEndomaxx_T3                       0.51542821
## SPAD.avg                                                0.71055300
## Biostimulant_TimePointSumagrow_T2                       0.51542821
## Biostimulant_TimePointPathway_T2                        0.51542821
## BiostimulantSumagrow                                    0.45335192
## Biostimulant_TimePointInocucor_T2                       0.51542821
## Biostimulant_TimePointControl_T2                        1.00000000
## BiostimulantPathway                                     0.45335192
## BiostimulantInocucor                                    0.45335192
## Biostimulant_TimePointEndomaxx_T2                       0.51542821
## BiostimulantEndomaxx                                    0.45335192
## Biostimulant_TimePointInocucor_T1                       0.51542821
## Biostimulant_TimePointSumagrow_T1                       0.06255752
## Biostimulant_TimePointPathway_T1                        0.12724239
## Biostimulant_TimePointControl_T1                        0.70280105
## Biostimulant_TimePointEndomaxx_T1                       0.86395099
## Biostimulant_TimePointPathway_T0                        0.92310273
## Biostimulant_TimePointSumagrow_T0                       0.89266416
## Biostimulant_TimePointInocucor_T0                       0.06255752
## Biostimulant_TimePointEndomaxx_T0                       0.93718225
##                                   BiostimulantPathway BiostimulantInocucor
## Leaves                                    0.211222326          0.231592502
## Tot                                       0.225802627          0.218743512
## Root                                      0.020993415          0.197864978
## Stems                                     0.233099163          0.007143058
## LDA1_biomass_vs_treatment                 0.592183861          0.549022196
## Area                                      0.545793957          0.596014664
## LDA1_area_vs_treatment                    0.545793957          0.596014664
## Biostimulant_TimePointSumagrow_T3         0.453351923          0.453351923
## Biostimulant_TimePointControl_T3          0.453351923          0.453351923
## Biostimulant_TimePointInocucor_T3         0.453351923          1.026891258
## Biostimulant_TimePointPathway_T3          1.026891258          0.453351923
## Biostimulant_TimePointEndomaxx_T3         0.453351923          0.453351923
## SPAD.avg                                  0.639186188          0.518819309
## Biostimulant_TimePointSumagrow_T2         0.453351923          0.453351923
## Biostimulant_TimePointPathway_T2          1.026891258          0.453351923
## BiostimulantSumagrow                      0.318059790          0.318059790
## Biostimulant_TimePointInocucor_T2         0.453351923          1.026891258
## Biostimulant_TimePointControl_T2          0.453351923          0.453351923
## BiostimulantPathway                       1.000000000          0.318059790
## BiostimulantInocucor                      0.318059790          1.000000000
## Biostimulant_TimePointEndomaxx_T2         0.453351923          0.453351923
## BiostimulantEndomaxx                      0.318059790          0.318059790
## Biostimulant_TimePointInocucor_T1         0.453351923          1.026891258
## Biostimulant_TimePointSumagrow_T1         0.005733727         -0.114524032
## Biostimulant_TimePointPathway_T1          1.175358918         -0.111853048
## Biostimulant_TimePointControl_T1         -0.114584388         -0.105534803
## Biostimulant_TimePointEndomaxx_T1        -0.114696589         -0.078477527
## Biostimulant_TimePointPathway_T0          0.458858014         -0.087571573
## Biostimulant_TimePointSumagrow_T0        -0.114680397         -0.093943654
## Biostimulant_TimePointInocucor_T0         0.005733727          0.459015303
## Biostimulant_TimePointEndomaxx_T0        -0.114684520         -0.082079474
##                                   Biostimulant_TimePointEndomaxx_T2
## Leaves                                                  0.194268115
## Tot                                                     0.182758507
## Root                                                    0.164054286
## Stems                                                  -0.028229295
## LDA1_biomass_vs_treatment                               0.515919081
## Area                                                    0.695747215
## LDA1_area_vs_treatment                                  0.695747215
## Biostimulant_TimePointSumagrow_T3                       0.515428211
## Biostimulant_TimePointControl_T3                        0.515428211
## Biostimulant_TimePointInocucor_T3                       0.515428211
## Biostimulant_TimePointPathway_T3                        0.515428211
## Biostimulant_TimePointEndomaxx_T3                       0.515428211
## SPAD.avg                                                0.514498402
## Biostimulant_TimePointSumagrow_T2                       0.515428211
## Biostimulant_TimePointPathway_T2                        0.515428211
## BiostimulantSumagrow                                    0.453351923
## Biostimulant_TimePointInocucor_T2                       0.515428211
## Biostimulant_TimePointControl_T2                        0.515428211
## BiostimulantPathway                                     0.453351923
## BiostimulantInocucor                                    0.453351923
## Biostimulant_TimePointEndomaxx_T2                       1.000000000
## BiostimulantEndomaxx                                    1.026891258
## Biostimulant_TimePointInocucor_T1                       0.515428211
## Biostimulant_TimePointSumagrow_T1                       0.760285183
## Biostimulant_TimePointPathway_T1                        0.003878776
## Biostimulant_TimePointControl_T1                        0.132008612
## Biostimulant_TimePointEndomaxx_T1                       0.057975685
## Biostimulant_TimePointPathway_T0                        0.096080454
## Biostimulant_TimePointSumagrow_T0                       0.086057018
## Biostimulant_TimePointInocucor_T0                       0.760285183
## Biostimulant_TimePointEndomaxx_T0                       0.093225514
##                                   BiostimulantEndomaxx
## Leaves                                      0.16011999
## Tot                                         0.17525013
## Root                                        0.20207966
## Stems                                      -0.03851243
## LDA1_biomass_vs_treatment                   0.56393470
## Area                                        0.54720190
## LDA1_area_vs_treatment                      0.54720190
## Biostimulant_TimePointSumagrow_T3           0.45335192
## Biostimulant_TimePointControl_T3            0.45335192
## Biostimulant_TimePointInocucor_T3           0.45335192
## Biostimulant_TimePointPathway_T3            0.45335192
## Biostimulant_TimePointEndomaxx_T3           1.02689126
## SPAD.avg                                    0.40553617
## Biostimulant_TimePointSumagrow_T2           0.45335192
## Biostimulant_TimePointPathway_T2            0.45335192
## BiostimulantSumagrow                        0.31805979
## Biostimulant_TimePointInocucor_T2           0.45335192
## Biostimulant_TimePointControl_T2            0.45335192
## BiostimulantPathway                         0.31805979
## BiostimulantInocucor                        0.31805979
## Biostimulant_TimePointEndomaxx_T2           1.02689126
## BiostimulantEndomaxx                        1.00000000
## Biostimulant_TimePointInocucor_T1           0.45335192
## Biostimulant_TimePointSumagrow_T1          -0.03622074
## Biostimulant_TimePointPathway_T1           -0.09243608
## Biostimulant_TimePointControl_T1            0.60311508
## Biostimulant_TimePointEndomaxx_T1           1.23362256
## Biostimulant_TimePointPathway_T0            0.76480168
## Biostimulant_TimePointSumagrow_T0           0.85753219
## Biostimulant_TimePointInocucor_T0          -0.03622074
## Biostimulant_TimePointEndomaxx_T0           1.27842588
##                                   Biostimulant_TimePointInocucor_T1
## Leaves                                                   0.09170286
## Tot                                                      0.07834586
## Root                                                     0.08324646
## Stems                                                   -0.14461099
## LDA1_biomass_vs_treatment                                0.43846561
## Area                                                     0.63552361
## LDA1_area_vs_treatment                                   0.63552361
## Biostimulant_TimePointSumagrow_T3                        0.51542821
## Biostimulant_TimePointControl_T3                         0.51542821
## Biostimulant_TimePointInocucor_T3                        0.51542821
## Biostimulant_TimePointPathway_T3                         0.51542821
## Biostimulant_TimePointEndomaxx_T3                        0.51542821
## SPAD.avg                                                 0.75444427
## Biostimulant_TimePointSumagrow_T2                        0.51542821
## Biostimulant_TimePointPathway_T2                         0.51542821
## BiostimulantSumagrow                                     0.45335192
## Biostimulant_TimePointInocucor_T2                        0.51542821
## Biostimulant_TimePointControl_T2                         0.51542821
## BiostimulantPathway                                      0.45335192
## BiostimulantInocucor                                     1.02689126
## Biostimulant_TimePointEndomaxx_T2                        0.51542821
## BiostimulantEndomaxx                                     0.45335192
## Biostimulant_TimePointInocucor_T1                        1.00000000
## Biostimulant_TimePointSumagrow_T1                       -0.03505166
## Biostimulant_TimePointPathway_T1                         0.06292161
## Biostimulant_TimePointControl_T1                        -0.05250850
## Biostimulant_TimePointEndomaxx_T1                       -0.05261211
## Biostimulant_TimePointPathway_T0                        -0.05262428
## Biostimulant_TimePointSumagrow_T0                       -0.05260559
## Biostimulant_TimePointInocucor_T0                       -0.03505166
## Biostimulant_TimePointEndomaxx_T0                       -0.05262623
##                                   Biostimulant_TimePointSumagrow_T1
## Leaves                                                  0.699182697
## Tot                                                     0.625862955
## Root                                                    0.676622855
## Stems                                                   0.556747101
## LDA1_biomass_vs_treatment                               0.771902116
## Area                                                    0.083816638
## LDA1_area_vs_treatment                                  0.243999664
## Biostimulant_TimePointSumagrow_T3                      -0.052438928
## Biostimulant_TimePointControl_T3                        0.713707800
## Biostimulant_TimePointInocucor_T3                       0.111957723
## Biostimulant_TimePointPathway_T3                       -0.020237671
## Biostimulant_TimePointEndomaxx_T3                      -0.052291586
## SPAD.avg                                                0.638669744
## Biostimulant_TimePointSumagrow_T2                       0.035169636
## Biostimulant_TimePointPathway_T2                        0.051659546
## BiostimulantSumagrow                                    0.458976367
## Biostimulant_TimePointInocucor_T2                       0.843269216
## Biostimulant_TimePointControl_T2                        0.062557522
## BiostimulantPathway                                     0.005733727
## BiostimulantInocucor                                   -0.114524032
## Biostimulant_TimePointEndomaxx_T2                       0.760285183
## BiostimulantEndomaxx                                   -0.036220738
## Biostimulant_TimePointInocucor_T1                      -0.035051664
## Biostimulant_TimePointSumagrow_T1                       1.000000000
## Biostimulant_TimePointPathway_T1                       -0.052631579
## Biostimulant_TimePointControl_T1                       -0.052631579
## Biostimulant_TimePointEndomaxx_T1                      -0.052631579
## Biostimulant_TimePointPathway_T0                       -0.052631579
## Biostimulant_TimePointSumagrow_T0                      -0.052631579
## Biostimulant_TimePointInocucor_T0                      -0.052631579
## Biostimulant_TimePointEndomaxx_T0                      -0.052631579
##                                   Biostimulant_TimePointPathway_T1
## Leaves                                                -0.060371369
## Tot                                                    0.450450390
## Root                                                   0.304314809
## Stems                                                 -0.022550517
## LDA1_biomass_vs_treatment                             -0.106270215
## Area                                                   0.132605457
## LDA1_area_vs_treatment                                 0.302496746
## Biostimulant_TimePointSumagrow_T3                     -0.052614282
## Biostimulant_TimePointControl_T3                       0.018430601
## Biostimulant_TimePointInocucor_T3                      0.508601066
## Biostimulant_TimePointPathway_T3                       0.702584843
## Biostimulant_TimePointEndomaxx_T3                     -0.052205736
## SPAD.avg                                               0.233784243
## Biostimulant_TimePointSumagrow_T2                      0.695264952
## Biostimulant_TimePointPathway_T2                       0.858919594
## BiostimulantSumagrow                                  -0.016093440
## Biostimulant_TimePointInocucor_T2                     -0.017230150
## Biostimulant_TimePointControl_T2                       0.127242394
## BiostimulantPathway                                    1.175358918
## BiostimulantInocucor                                  -0.111853048
## Biostimulant_TimePointEndomaxx_T2                      0.003878776
## BiostimulantEndomaxx                                  -0.092436078
## Biostimulant_TimePointInocucor_T1                      0.062921609
## Biostimulant_TimePointSumagrow_T1                     -0.052631579
## Biostimulant_TimePointPathway_T1                       1.000000000
## Biostimulant_TimePointControl_T1                      -0.037879867
## Biostimulant_TimePointEndomaxx_T1                      0.185867693
## Biostimulant_TimePointPathway_T0                       0.073345063
## Biostimulant_TimePointSumagrow_T0                      0.037071373
## Biostimulant_TimePointInocucor_T0                     -0.052631579
## Biostimulant_TimePointEndomaxx_T0                      0.113607469
##                                   Biostimulant_TimePointControl_T1
## Leaves                                                  0.44454594
## Tot                                                     0.68610099
## Root                                                    0.84565505
## Stems                                                   0.75548059
## LDA1_biomass_vs_treatment                               0.01715083
## Area                                                    0.88853235
## LDA1_area_vs_treatment                                  0.01096427
## Biostimulant_TimePointSumagrow_T3                      -0.04379979
## Biostimulant_TimePointControl_T3                        0.06427478
## Biostimulant_TimePointInocucor_T3                       0.79321433
## Biostimulant_TimePointPathway_T3                       -0.05190465
## Biostimulant_TimePointEndomaxx_T3                      -0.04365947
## SPAD.avg                                                0.38956393
## Biostimulant_TimePointSumagrow_T2                       0.61044040
## Biostimulant_TimePointPathway_T2                       -0.05255685
## BiostimulantSumagrow                                   -0.10566186
## Biostimulant_TimePointInocucor_T2                       0.08445510
## Biostimulant_TimePointControl_T2                        0.70280105
## BiostimulantPathway                                    -0.11458439
## BiostimulantInocucor                                   -0.10553480
## Biostimulant_TimePointEndomaxx_T2                       0.13200861
## BiostimulantEndomaxx                                    0.60311508
## Biostimulant_TimePointInocucor_T1                      -0.05250850
## Biostimulant_TimePointSumagrow_T1                      -0.05263158
## Biostimulant_TimePointPathway_T1                       -0.03787987
## Biostimulant_TimePointControl_T1                        1.00000000
## Biostimulant_TimePointEndomaxx_T1                      -0.05263158
## Biostimulant_TimePointPathway_T0                       -0.05263158
## Biostimulant_TimePointSumagrow_T0                      -0.05263158
## Biostimulant_TimePointInocucor_T0                      -0.05263158
## Biostimulant_TimePointEndomaxx_T0                      -0.05263158
##                                   Biostimulant_TimePointEndomaxx_T1
## Leaves                                                   0.53233412
## Tot                                                      0.78757103
## Root                                                     0.68139120
## Stems                                                    0.77009775
## LDA1_biomass_vs_treatment                               -0.05805714
## Area                                                     0.58601563
## LDA1_area_vs_treatment                                  -0.02740427
## Biostimulant_TimePointSumagrow_T3                       -0.01792098
## Biostimulant_TimePointControl_T3                         0.02945594
## Biostimulant_TimePointInocucor_T3                        0.70692574
## Biostimulant_TimePointPathway_T3                        -0.05205926
## Biostimulant_TimePointEndomaxx_T3                       -0.01724260
## SPAD.avg                                                 0.28043536
## Biostimulant_TimePointSumagrow_T2                        0.76283252
## Biostimulant_TimePointPathway_T2                        -0.05260214
## BiostimulantSumagrow                                    -0.07906615
## Biostimulant_TimePointInocucor_T2                        0.05172942
## Biostimulant_TimePointControl_T2                         0.86395099
## BiostimulantPathway                                     -0.11469659
## BiostimulantInocucor                                    -0.07847753
## Biostimulant_TimePointEndomaxx_T2                        0.05797568
## BiostimulantEndomaxx                                     1.23362256
## Biostimulant_TimePointInocucor_T1                       -0.05261211
## Biostimulant_TimePointSumagrow_T1                       -0.05263158
## Biostimulant_TimePointPathway_T1                         0.18586769
## Biostimulant_TimePointControl_T1                        -0.05263158
## Biostimulant_TimePointEndomaxx_T1                        1.00000000
## Biostimulant_TimePointPathway_T0                        -0.05263158
## Biostimulant_TimePointSumagrow_T0                       -0.05263158
## Biostimulant_TimePointInocucor_T0                       -0.05263158
## Biostimulant_TimePointEndomaxx_T0                       -0.05263158
##                                   Biostimulant_TimePointPathway_T0
## Leaves                                                 0.675264423
## Tot                                                    0.665134628
## Root                                                   0.618048832
## Stems                                                  0.583524966
## LDA1_biomass_vs_treatment                             -0.083451608
## Area                                                   0.044837995
## LDA1_area_vs_treatment                                -0.340011255
## Biostimulant_TimePointSumagrow_T3                     -0.026526753
## Biostimulant_TimePointControl_T3                       0.076284710
## Biostimulant_TimePointInocucor_T3                      0.548885265
## Biostimulant_TimePointPathway_T3                      -0.052602952
## Biostimulant_TimePointEndomaxx_T3                     -0.025868689
## SPAD.avg                                               0.008915816
## Biostimulant_TimePointSumagrow_T2                      0.403323976
## Biostimulant_TimePointPathway_T2                      -0.052605175
## BiostimulantSumagrow                                  -0.087712117
## Biostimulant_TimePointInocucor_T2                      0.086489902
## Biostimulant_TimePointControl_T2                       0.923102735
## BiostimulantPathway                                    0.458858014
## BiostimulantInocucor                                  -0.087571573
## Biostimulant_TimePointEndomaxx_T2                      0.096080454
## BiostimulantEndomaxx                                   0.764801678
## Biostimulant_TimePointInocucor_T1                     -0.052624276
## Biostimulant_TimePointSumagrow_T1                     -0.052631579
## Biostimulant_TimePointPathway_T1                       0.073345063
## Biostimulant_TimePointControl_T1                      -0.052631579
## Biostimulant_TimePointEndomaxx_T1                     -0.052631579
## Biostimulant_TimePointPathway_T0                       1.000000000
## Biostimulant_TimePointSumagrow_T0                     -0.052631579
## Biostimulant_TimePointInocucor_T0                     -0.052631579
## Biostimulant_TimePointEndomaxx_T0                     -0.052631579
##                                   Biostimulant_TimePointSumagrow_T0
## Leaves                                                   0.45171076
## Tot                                                      0.73256934
## Root                                                     0.68105724
## Stems                                                    0.69690528
## LDA1_biomass_vs_treatment                               -0.08708204
## Area                                                     0.38003179
## LDA1_area_vs_treatment                                  -0.33501550
## Biostimulant_TimePointSumagrow_T3                       -0.03270557
## Biostimulant_TimePointControl_T3                         0.04773485
## Biostimulant_TimePointInocucor_T3                        0.69691538
## Biostimulant_TimePointPathway_T3                        -0.05233390
## Biostimulant_TimePointEndomaxx_T3                       -0.03230415
## SPAD.avg                                                -0.12508802
## Biostimulant_TimePointSumagrow_T2                        0.85038642
## Biostimulant_TimePointPathway_T2                        -0.05259743
## BiostimulantSumagrow                                     0.47932140
## Biostimulant_TimePointInocucor_T2                        0.06700611
## Biostimulant_TimePointControl_T2                         0.89266416
## BiostimulantPathway                                     -0.11468040
## BiostimulantInocucor                                    -0.09394365
## Biostimulant_TimePointEndomaxx_T2                        0.08605702
## BiostimulantEndomaxx                                     0.85753219
## Biostimulant_TimePointInocucor_T1                       -0.05260559
## Biostimulant_TimePointSumagrow_T1                       -0.05263158
## Biostimulant_TimePointPathway_T1                         0.03707137
## Biostimulant_TimePointControl_T1                        -0.05263158
## Biostimulant_TimePointEndomaxx_T1                       -0.05263158
## Biostimulant_TimePointPathway_T0                        -0.05263158
## Biostimulant_TimePointSumagrow_T0                        1.00000000
## Biostimulant_TimePointInocucor_T0                       -0.05263158
## Biostimulant_TimePointEndomaxx_T0                       -0.05263158
##                                   Biostimulant_TimePointInocucor_T0
## Leaves                                                  0.583503928
## Tot                                                     0.550206651
## Root                                                    0.607035263
## Stems                                                   0.510865654
## LDA1_biomass_vs_treatment                               0.709799496
## Area                                                   -0.259701868
## LDA1_area_vs_treatment                                 -0.099518842
## Biostimulant_TimePointSumagrow_T3                      -0.052438928
## Biostimulant_TimePointControl_T3                        0.713707800
## Biostimulant_TimePointInocucor_T3                       0.111957723
## Biostimulant_TimePointPathway_T3                       -0.020237671
## Biostimulant_TimePointEndomaxx_T3                      -0.052291586
## SPAD.avg                                                0.147441171
## Biostimulant_TimePointSumagrow_T2                       0.035169636
## Biostimulant_TimePointPathway_T2                        0.051659546
## BiostimulantSumagrow                                   -0.114562967
## Biostimulant_TimePointInocucor_T2                       0.843269216
## Biostimulant_TimePointControl_T2                        0.062557522
## BiostimulantPathway                                     0.005733727
## BiostimulantInocucor                                    0.459015303
## Biostimulant_TimePointEndomaxx_T2                       0.760285183
## BiostimulantEndomaxx                                   -0.036220738
## Biostimulant_TimePointInocucor_T1                      -0.035051664
## Biostimulant_TimePointSumagrow_T1                      -0.052631579
## Biostimulant_TimePointPathway_T1                       -0.052631579
## Biostimulant_TimePointControl_T1                       -0.052631579
## Biostimulant_TimePointEndomaxx_T1                      -0.052631579
## Biostimulant_TimePointPathway_T0                       -0.052631579
## Biostimulant_TimePointSumagrow_T0                      -0.052631579
## Biostimulant_TimePointInocucor_T0                       1.000000000
## Biostimulant_TimePointEndomaxx_T0                      -0.052631579
##                                   Biostimulant_TimePointEndomaxx_T0
## Leaves                                                   0.72625469
## Tot                                                      0.62420587
## Root                                                     0.59059149
## Stems                                                    0.56503191
## LDA1_biomass_vs_treatment                               -0.08952663
## Area                                                    -0.02269113
## LDA1_area_vs_treatment                                  -0.33401142
## Biostimulant_TimePointSumagrow_T3                       -0.02123148
## Biostimulant_TimePointControl_T3                         0.07993096
## Biostimulant_TimePointInocucor_T3                        0.51904818
## Biostimulant_TimePointPathway_T3                        -0.05261217
## Biostimulant_TimePointEndomaxx_T3                       -0.02041589
## SPAD.avg                                                -0.22301769
## Biostimulant_TimePointSumagrow_T2                        0.30932239
## Biostimulant_TimePointPathway_T2                        -0.05260617
## BiostimulantSumagrow                                    -0.08221434
## Biostimulant_TimePointInocucor_T2                        0.08880906
## Biostimulant_TimePointControl_T2                         0.93718225
## BiostimulantPathway                                     -0.11468452
## BiostimulantInocucor                                    -0.08207947
## Biostimulant_TimePointEndomaxx_T2                        0.09322551
## BiostimulantEndomaxx                                     1.27842588
## Biostimulant_TimePointInocucor_T1                       -0.05262623
## Biostimulant_TimePointSumagrow_T1                       -0.05263158
## Biostimulant_TimePointPathway_T1                         0.11360747
## Biostimulant_TimePointControl_T1                        -0.05263158
## Biostimulant_TimePointEndomaxx_T1                       -0.05263158
## Biostimulant_TimePointPathway_T0                        -0.05263158
## Biostimulant_TimePointSumagrow_T0                       -0.05263158
## Biostimulant_TimePointInocucor_T0                       -0.05263158
## Biostimulant_TimePointEndomaxx_T0                        1.00000000

Tukey’s HSD

library(agricolae)

formula <- lda_variables %>%
  dplyr::select(-X.SampleID, -Biostimulant, -TimePoint, -Rep, -Biostimulant_TimePoint) %>%
  colnames(.) %>%
  paste(collapse = '+') %>%
  paste0("~Biostimulant*TimePoint")

model <- aov(formula = as.formula(formula),data = lda_variables)
HSD.test(model, c('Biostimulant','TimePoint'), group=TRUE, console = TRUE)
## 
## Study: model ~ c("Biostimulant", "TimePoint")
## 
## HSD Test for Area + SPAD.avg + Leaves + Stems + Root + Tot + LDA1_area_vs_treatment + LDA1_biomass_vs_treatment 
## 
## Mean Square Error:  59860423 
## 
## Biostimulant:TimePoint,  means
## 
##             Area...SPAD.avg...Leaves...Stems...Root...Tot...LDA1_area_vs_treatment...LDA1_biomass_vs_treatment
## Control:T0                                                                                            1423.980
## Control:T1                                                                                           11221.281
## Control:T2                                                                                           27919.769
## Control:T3                                                                                           90523.541
## Endomaxx:T0                                                                                           1233.455
## Endomaxx:T1                                                                                          10655.031
## Endomaxx:T2                                                                                          28458.386
## Endomaxx:T3                                                                                          79857.140
## Inocucor:T0                                                                                           1351.726
## Inocucor:T1                                                                                          13246.916
## Inocucor:T2                                                                                          31142.004
## Inocucor:T3                                                                                          86171.488
## Pathway:T0                                                                                            1366.712
## Pathway:T1                                                                                           11859.052
## Pathway:T2                                                                                           33673.238
## Pathway:T3                                                                                           86753.210
## Sumagrow:T0                                                                                           1435.490
## Sumagrow:T1                                                                                          12869.290
## Sumagrow:T2                                                                                          32320.106
## Sumagrow:T3                                                                                          87026.510
##                    std r        Min        Max
## Control:T0    209.6647 6  1127.5490   1753.384
## Control:T1   1608.8545 6  9160.6688  13159.208
## Control:T2  11827.2449 6  6182.1711  38148.368
## Control:T3  11282.4934 6 71122.3143 103970.407
## Endomaxx:T0   336.6803 6   749.9229   1716.401
## Endomaxx:T1  2569.9170 6  6252.4763  13009.854
## Endomaxx:T2  8004.2091 6 12648.9982  35096.270
## Endomaxx:T3  6390.8685 6 68436.2083  87294.333
## Inocucor:T0   164.4797 6  1149.8716   1536.325
## Inocucor:T1   757.6288 6 12368.0014  14399.580
## Inocucor:T2  3135.2540 6 27487.8944  35814.160
## Inocucor:T3 22391.1987 6 64551.6091 126706.771
## Pathway:T0    100.4050 6  1244.5531   1510.122
## Pathway:T1   2537.1205 6  8119.7624  15480.086
## Pathway:T2   7134.4230 6 26832.5593  46420.858
## Pathway:T3  10122.5326 6 74033.4369 104154.662
## Sumagrow:T0   312.1367 6   892.4997   1730.269
## Sumagrow:T1   731.6217 6 11765.6893  13642.374
## Sumagrow:T2  8501.1937 6 22234.4478  42830.223
## Sumagrow:T3  8441.3829 6 75258.8619  98553.018
## 
## Alpha: 0.05 ; DF Error: 100 
## Critical Value of Studentized Range: 5.148913 
## 
## Minimun Significant Difference: 16263.34 
## 
## Treatments with the same letter are not significantly different.
## 
##             Area + SPAD.avg + Leaves + Stems + Root + Tot + LDA1_area_vs_treatment + LDA1_biomass_vs_treatment
## Control:T3                                                                                           90523.541
## Sumagrow:T3                                                                                          87026.510
## Pathway:T3                                                                                           86753.210
## Inocucor:T3                                                                                          86171.488
## Endomaxx:T3                                                                                          79857.140
## Pathway:T2                                                                                           33673.238
## Sumagrow:T2                                                                                          32320.106
## Inocucor:T2                                                                                          31142.004
## Endomaxx:T2                                                                                          28458.386
## Control:T2                                                                                           27919.769
## Inocucor:T1                                                                                          13246.916
## Sumagrow:T1                                                                                          12869.290
## Pathway:T1                                                                                           11859.052
## Control:T1                                                                                           11221.281
## Endomaxx:T1                                                                                          10655.031
## Sumagrow:T0                                                                                           1435.490
## Control:T0                                                                                            1423.980
## Pathway:T0                                                                                            1366.712
## Inocucor:T0                                                                                           1351.726
## Endomaxx:T0                                                                                           1233.455
##             groups
## Control:T3       a
## Sumagrow:T3      a
## Pathway:T3       a
## Inocucor:T3      a
## Endomaxx:T3      a
## Pathway:T2       b
## Sumagrow:T2      b
## Inocucor:T2      b
## Endomaxx:T2     bc
## Control:T2     bcd
## Inocucor:T1    cde
## Sumagrow:T1    cde
## Pathway:T1      de
## Control:T1       e
## Endomaxx:T1      e
## Sumagrow:T0      e
## Control:T0       e
## Pathway:T0       e
## Inocucor:T0      e
## Endomaxx:T0      e
readr::write_tsv(lda_variables,'mappings_LDA.txt')

Phyloseq

Import files

library(phyloseq)
library(vegan)
## Loading required package: permute
## Loading required package: lattice
## This is vegan 2.5-2
library(ape)
## 
## Attaching package: 'ape'
## The following object is masked from 'package:agricolae':
## 
##     consensus
library(dummies)
## dummies-1.5.6 provided by Decision Patterns
#File Paths
biom_path <- file.path('merged/biom/table_wo_chl_mit.biom')
tree_path <- file.path('merged/biom/tree.nwk')
DESEq2_path <- file.path('merged/biom/DESeq2_w_tax.biom')
map_path <- file.path('mappings_LDA.txt')

#Import to phyloseq table and merge into phyloseq objects
table <- import_biom(BIOMfilename = biom_path,
                      #refseqfilename = repseqfile,
                      parseFunction = parse_taxonomy_default, 
                      parallel = T)
## Warning in strsplit(conditionMessage(e), "\n"): input string 1 is invalid
## in this locale
tax_table(table) <-tax_table(table)[,1:7]
DESEq2.table <- import_biom(BIOMfilename = DESEq2_path,
                      #refseqfilename = repseqfile,
                      parseFunction = parse_taxonomy_default, 
                      parallel = T)
## Warning in strsplit(conditionMessage(e), "\n"): input string 1 is invalid
## in this locale
tax_table(DESEq2.table) <-tax_table(DESEq2.table)[,1:7]
metadata <- import_qiime_sample_data(map_path)
tree <- read_tree(tree_path)


phylobj <- merge_phyloseq(table, metadata, tree)
DESEq2.phylobj <- merge_phyloseq(DESEq2.table, metadata, tree)

#Adjust taxonomy names (to harmonize betwen UNITE and SILVA databases)
tax_table(phylobj) <- gsub(".*__", "", tax_table(phylobj))
colnames(tax_table(phylobj)) <- c("Kingdom", "Phylum", "Class", 
                    "Order", "Family", "Genus", "Species")

tax_table(DESEq2.phylobj) <- gsub(".*__", "", tax_table(DESEq2.phylobj))
colnames(tax_table(DESEq2.phylobj)) <- c("Kingdom", "Phylum", "Class", 
                    "Order", "Family", "Genus", "Species")

#Log-transform sample counts (if needed)
pslog <- transform_sample_counts(phylobj,function(x){log(1 + x)})

Phylum barplots

phyloglom <- tax_glom(phylobj, 'Phylum')

phylorel <- transform_sample_counts(phyloglom, function(x){100*x/sum(x)})

phylorel %>%
  filter_taxa(function(x){ mean(x) > 1}, TRUE) %>%
  psmelt()%>%
  group_by(Biostimulant, TimePoint, Phylum) %>%
  summarize(mean=mean(Abundance)) %>%
  ggplot() +
  aes(x=Biostimulant, y=mean, fill=Phylum, color=Phylum) +
  geom_bar(stat='identity') +
  facet_grid(.~TimePoint) +
  ylab('Relative abundance') +
  theme_igray() +
  theme(axis.text.x  = element_text(angle=90, vjust=0.5)) +
  scale_fill_pander() + scale_color_pander()

Alpha Diveristy

plot_richness(phylobj, 
              x = "Biostimulant", 
              color = "TimePoint", 
              measures = c('Observed', 'Chao1', 'Shannon')) + 
  geom_boxplot(alpha=.9) + theme_igray() +
  theme(axis.text.x  = element_text(angle=90, vjust=0.5))+
  scale_color_pander() + scale_fill_pander()
## Warning: Removed 240 rows containing missing values (geom_errorbar).

##Rarefaction plots

source('https://raw.githubusercontent.com/mahendra-mariadassou/phyloseq-extended/master/R/richness.R')
## Loading required package: parallel
ggrare(phylobj, step = 100, color = "Biostimulant", label = "Sample", se = FALSE) + 
  facet_wrap(~TimePoint) + guides(label=FALSE) + theme_igray() + scale_color_pander()
## rarefying sample B2CT0-1
## rarefying sample B2CT0-2
## rarefying sample B2CT0-3
## rarefying sample B2CT0-4
## rarefying sample B2CT0-5
## rarefying sample B2CT0-6
## rarefying sample B2CT1-1
## rarefying sample B2CT1-2
## rarefying sample B2CT1-3
## rarefying sample B2CT1-4
## rarefying sample B2CT1-5
## rarefying sample B2CT1-6
## rarefying sample B2CT2-1
## rarefying sample B2CT2-2
## rarefying sample B2CT2-3
## rarefying sample B2CT2-4
## rarefying sample B2CT2-5
## rarefying sample B2CT2-6
## rarefying sample B2CT3-1
## rarefying sample B2CT3-2
## rarefying sample B2CT3-3
## rarefying sample B2CT3-4
## rarefying sample B2CT3-5
## rarefying sample B2CT3-6
## rarefying sample B2ET0-1
## rarefying sample B2ET0-2
## rarefying sample B2ET0-3
## rarefying sample B2ET0-4
## rarefying sample B2ET0-5
## rarefying sample B2ET0-6
## rarefying sample B2ET1-1
## rarefying sample B2ET1-2
## rarefying sample B2ET1-3
## rarefying sample B2ET1-4
## rarefying sample B2ET1-5
## rarefying sample B2ET1-6
## rarefying sample B2ET2-1
## rarefying sample B2ET2-2
## rarefying sample B2ET2-3
## rarefying sample B2ET2-4
## rarefying sample B2ET2-5
## rarefying sample B2ET2-6
## rarefying sample B2ET3-1
## rarefying sample B2ET3-2
## rarefying sample B2ET3-3
## rarefying sample B2ET3-4
## rarefying sample B2ET3-5
## rarefying sample B2ET3-6
## rarefying sample B2IT0-1
## rarefying sample B2IT0-2
## rarefying sample B2IT0-3
## rarefying sample B2IT0-4
## rarefying sample B2IT0-5
## rarefying sample B2IT0-6
## rarefying sample B2IT1-1
## rarefying sample B2IT1-2
## rarefying sample B2IT1-3
## rarefying sample B2IT1-4
## rarefying sample B2IT1-5
## rarefying sample B2IT1-6
## rarefying sample B2IT2-1
## rarefying sample B2IT2-2
## rarefying sample B2IT2-3
## rarefying sample B2IT2-4
## rarefying sample B2IT2-5
## rarefying sample B2IT2-6
## rarefying sample B2IT3-1
## rarefying sample B2IT3-2
## rarefying sample B2IT3-3
## rarefying sample B2IT3-4
## rarefying sample B2IT3-5
## rarefying sample B2IT3-6
## rarefying sample B2PT0-1
## rarefying sample B2PT0-2
## rarefying sample B2PT0-3
## rarefying sample B2PT0-4
## rarefying sample B2PT0-5
## rarefying sample B2PT0-6
## rarefying sample B2PT1-1
## rarefying sample B2PT1-2
## rarefying sample B2PT1-3
## rarefying sample B2PT1-4
## rarefying sample B2PT1-5
## rarefying sample B2PT1-6
## rarefying sample B2PT2-1
## rarefying sample B2PT2-2
## rarefying sample B2PT2-3
## rarefying sample B2PT2-4
## rarefying sample B2PT2-5
## rarefying sample B2PT2-6
## rarefying sample B2PT3-1
## rarefying sample B2PT3-2
## rarefying sample B2PT3-3
## rarefying sample B2PT3-4
## rarefying sample B2PT3-5
## rarefying sample B2PT3-6
## rarefying sample B2ST0-1
## rarefying sample B2ST0-2
## rarefying sample B2ST0-3
## rarefying sample B2ST0-4
## rarefying sample B2ST0-5
## rarefying sample B2ST0-6
## rarefying sample B2ST1-1
## rarefying sample B2ST1-2
## rarefying sample B2ST1-3
## rarefying sample B2ST1-4
## rarefying sample B2ST1-5
## rarefying sample B2ST1-6
## rarefying sample B2ST2-1
## rarefying sample B2ST2-2
## rarefying sample B2ST2-3
## rarefying sample B2ST2-4
## rarefying sample B2ST2-5
## rarefying sample B2ST2-6
## rarefying sample B2ST3-1
## rarefying sample B2ST3-2
## rarefying sample B2ST3-3
## rarefying sample B2ST3-4
## rarefying sample B2ST3-5
## rarefying sample B2ST3-6

PCoA with Bray distances

PCoA.bray.ord<- ordinate(DESEq2.phylobj, "PCoA", distance = "bray")
PCoA.bray.plot <- plot_ordination(physeq = pslog,
                                        ordination = PCoA.bray.ord,
                                        color = "Biostimulant",
                                        shape = "TimePoint",
                                        title = "PCoA (Bray-Curtis distances)")
PCoA.bray.plot +geom_point(size=3) + scale_color_pander()

##PCoA with Unifrac distances

PCoA.bray.ord<- ordinate(pslog, "PCoA", distance = "unifrac")
PCoA.bray.plot <- plot_ordination(physeq = pslog,
                                        ordination = PCoA.bray.ord,
                                        color = "Biostimulant",
                                        shape = "TimePoint",
                                        title = "PCoA (Unifrac distances)")
PCoA.bray.plot +geom_point(size=3) + scale_color_pander()

If you want to run PCoA with ALL possible distances (adapted from Phyloseq website https://joey711.github.io/phyloseq/distance.html) (requires a lot of memory)

dist_methods <- unlist(distanceMethodList)
dist_methods = dist_methods[-which(dist_methods=="ANY")]

plist <- vector("list", length(dist_methods))
names(plist) = dist_methods
for( i in dist_methods ){
    # Calculate distance matrix
    iDist <- distance(pslog, method=i)
    # Calculate ordination
    iMDS  <- ordinate(pslog, "PCoA", distance=iDist)
    ## Make plot
    # Don't carry over previous plot (if error, p will be blank)
    p <- NULL
    # Create plot, store as temp variable, p
    p <- plot_ordination(pslog, iMDS, color="Biostimulant", shape="TimePoint")
    # Add title to each plot
    p <- p + ggtitle(paste("PCoA using distance method ", i, sep=""))
    # Save the graphic to file.
    plist[[i]] = p
}

df = ldply(plist, function(x) x$data)
names(df)[1] <- "distance"
p = ggplot(df, aes(Axis.1, Axis.2, color=Biostimulant, shape=TimePoint))
p = p + geom_point(size=3, alpha=0.5)
p = p + facet_wrap(~distance, scales="free")
p = p + ggtitle("PCoA on various distance metrics")
p

NMDS

NMDS_bray_b1 <- ordinate(DESEq2.phylobj, method = "NMDS", distance = "bray")
## Wisconsin double standardization
## Run 0 stress 0.03550533 
## Run 1 stress 0.04837901 
## Run 2 stress 0.03673057 
## Run 3 stress 0.03601908 
## Run 4 stress 0.03610428 
## Run 5 stress 0.05136521 
## Run 6 stress 0.03624341 
## Run 7 stress 0.03911558 
## Run 8 stress 0.04939574 
## Run 9 stress 0.03635159 
## Run 10 stress 0.04983131 
## Run 11 stress 0.0365035 
## Run 12 stress 0.04808414 
## Run 13 stress 0.03634461 
## Run 14 stress 0.03734913 
## Run 15 stress 0.05002115 
## Run 16 stress 0.05185519 
## Run 17 stress 0.03556772 
## ... Procrustes: rmse 0.001972417  max resid 0.01388614 
## Run 18 stress 0.03580298 
## ... Procrustes: rmse 0.0062454  max resid 0.04616961 
## Run 19 stress 0.05155175 
## Run 20 stress 0.05777514 
## *** No convergence -- monoMDS stopping criteria:
##     17: no. of iterations >= maxit
##      3: stress ratio > sratmax
stressplot(NMDS_bray_b1)

plot_NMDS_bray_b1 <- plot_ordination(pslog,
                                  NMDS_bray_b1,
                                  type="samples",
                                  color="Biostimulant",
                                  shape="TimePoint")
plot_NMDS_bray_b1 + ggtitle("NMDS using Bray-Curtis distances") + scale_color_pander()

##PERMANOVA

DESEq2.df = as(sample_data(DESEq2.phylobj), "data.frame")
DESEq2.distbray = distance(DESEq2.phylobj, "bray")
adonis(DESEq2.distbray ~ TimePoint*Biostimulant, DESEq2.df)
## 
## Call:
## adonis(formula = DESEq2.distbray ~ TimePoint * Biostimulant,      data = DESEq2.df) 
## 
## Permutation: free
## Number of permutations: 999
## 
## Terms added sequentially (first to last)
## 
##                         Df SumsOfSqs   MeanSqs F.Model      R2 Pr(>F)    
## TimePoint                3  0.024838 0.0082795  5.7980 0.13002  0.001 ***
## Biostimulant             4  0.005316 0.0013289  0.9306 0.02783  0.653    
## TimePoint:Biostimulant  12  0.018086 0.0015072  1.0555 0.09467  0.215    
## Residuals              100  0.142798 0.0014280         0.74748           
## Total                  119  0.191038                   1.00000           
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

ANOSIM and beta dispersivity

attach(metadata)
NMDS_anosim = anosim(DESEq2.distbray, Biostimulant)
summary(NMDS_anosim)
## 
## Call:
## anosim(x = DESEq2.distbray, grouping = Biostimulant) 
## Dissimilarity: bray 
## 
## ANOSIM statistic R: 0.003099 
##       Significance: 0.254 
## 
## Permutation: free
## Number of permutations: 999
## 
## Upper quantiles of permutations (null model):
##     90%     95%   97.5%     99% 
## 0.00778 0.01135 0.01480 0.01963 
## 
## Dissimilarity ranks between and within classes:
##           0%     25%    50%     75% 100%    N
## Between    1 1795.50 3576.5 5360.50 7140 5760
## Control  114 2125.00 4360.0 5730.75 7139  276
## Endomaxx 771 2953.25 4398.5 6281.00 7135  276
## Inocucor   8 1840.50 3811.5 5355.00 7059  276
## Pathway   45 1690.75 3539.5 5610.25 7124  276
## Sumagrow  13  906.25 2010.0 3533.75 5491  276
plot(NMDS_anosim)

beta <- betadisper(DESEq2.distbray, DESEq2.df$Biostimulant) 
permutest(beta)
## 
## Permutation test for homogeneity of multivariate dispersions
## Permutation: free
## Number of permutations: 999
## 
## Response: Distances
##            Df   Sum Sq    Mean Sq      F N.Perm Pr(>F)  
## Groups      4 0.004390 0.00109760 2.6446    999  0.037 *
## Residuals 115 0.047728 0.00041503                       
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

CCA with VIF selection

source('https://andreanuzzo.github.io/Strausslab/vif.cca.bw_sel.R')

vifvariables = lda_variables 

vifvariables = vifvariables[,which(names(vifvariables)!='LinkerPrimerSequence' &
                        names(vifvariables)!='BarcodeSequence' &
                        names(vifvariables)!='Description' &
                        names(vifvariables)!='InputFileName' &
                        names(vifvariables)!='X.SampleID' &
                        names(vifvariables)!='Batch' &
                        names(vifvariables)!='Rep')]

cca_vif <- vif.cca.bw_sel(DESEq2.phylobj, vifvariables, threshold = 5)
## [1] "Removing variable Root with VIF:167.88"
## [1] "Removing variable LDA1_biomass_vs_treatment with VIF:130.9"
## [1] "Removing variable Stems with VIF:79.5"
## [1] "Removing variable TimePointT3 with VIF:25.48"
## [1] "Removing variable Leaves with VIF:22.74"
## [1] "Removing variable BiostimulantSumagrow with VIF:11.68"
plot(cca_vif)

anova(cca_vif)

CCA plot

cca_plot <- plot_ordination(physeq = pslog, 
                            ordination = cca_vif,
                            type= 'split',
                            color = "Biostimulant",
                            label = 'Phylum'
                            ) + aes(shape = TimePoint) + 
geom_point(aes(colour = Biostimulant))

cca_arrowmat <- scores(cca_vif, display = "bp") 
cca_arrowdf <- data.frame(labels = rownames(cca_arrowmat), cca_arrowmat)
cca_arrow_map <- aes(xend = CCA1, 
                 yend = CCA2, 
                 x = 0, 
                 y = 0, 
                 color = NULL,
                 shape = NULL)

cca_label_map <- aes(x = 1.3 * CCA1, 
                 y = 1.3 * CCA2, 
                 color = NULL, 
                 label = labels,
                 shape = NULL)

cca_arrowhead = arrow(length = unit(0.02, "npc"))
cca_plot + geom_segment(mapping = cca_arrow_map, size = .5, data = cca_arrowdf, 
               color = "black", arrow = cca_arrowhead) +
  geom_text(mapping = cca_label_map, size = 2, data = cca_arrowdf) +
  ggtitle("CCA with Bray-Curtis distances") + scale_color_pander()
## Warning: Removed 10215 rows containing missing values (geom_point).

##CAP Plot

source('https://andreanuzzo.github.io/Strausslab/vif.cap.bw_sel.R')
cap_vif <- vif.cap.bw_sel(pslog, vifvariables, threshold = 5)
## [1] "Removing variable Root with VIF:167.19"
## [1] "Removing variable LDA1_biomass_vs_treatment with VIF:130.4"
## [1] "Removing variable Stems with VIF:79.59"
## [1] "Removing variable TimePointT3 with VIF:25.45"
## [1] "Removing variable Leaves with VIF:22.6"
## [1] "Removing variable BiostimulantSumagrow with VIF:11.48"
cap_plot <- plot_ordination(physeq = pslog, 
                            ordination = cap_vif,
                            type= 'split',
                            color = "Biostimulant", 
                            label = 'Phylum',
                            axes = c(1,2)) + 
  aes(shape = TimePoint) +
  geom_point(aes(colour = Biostimulant)) + scale_color_pander()
arrowmat <- scores(cap_vif, display = "bp") 
arrowdf <- data.frame(labels = rownames(arrowmat), arrowmat)
arrow_map <- aes(xend = CAP1, 
                 yend = CAP2, 
                 x = 0, 
                 y = 0, 
                 color = NULL,
                 shape = NULL)
label_map <- aes(x = 1.3 * CAP1, 
                 y = 1.3 * CAP2, 
                 color = NULL, 
                 label = labels,
                 shape = NULL)
arrowhead = arrow(length = unit(0.02, "npc"))
cap_plot + 
  geom_segment(mapping = arrow_map, size = .5, data = arrowdf, 
               color = "#0C3D4C", arrow = arrowhead) +
  geom_text(mapping = label_map, size = 2, data = arrowdf, show.legend = FALSE) +
  ggtitle("Constrained Analysis with Bray-Curtis distances") 
## Warning: Removed 10215 rows containing missing values (geom_point).

Machine Learning in Qiime2

Input preparation

#!/bin/bash

# ----------------SLURM Parameters----------------

#SBATCH -J q2_model_input
#SBATCH --time=1:00:00
#SBATCH --ntasks=1
#SBATCH -D /ufrc/strauss/andrea.nuzzo/projects/Workshop/merged
#SBATCH -o logs/4_prepare_input_%j.out
#SBATCH -A strauss

date;hostname;pwd

################################################################################
#
# This script will prepare the input table for the modelling 
# 
################################################################################

# ----------------Housekeeping--------------------
rm -r models
mkdir models
mkdir models/input
cd models

# ----------------Load Modules--------------------
module load qiime2

# ------------------Commands----------------------
qiime taxa filter-table \
  --i-table ../features/table.qza \
  --i-taxonomy ../features/taxonomy.qza \
  --p-exclude chloroplast,mithocondria \
  --o-filtered-table input/inptab.qza

cd input

#Check the merged table

qiime feature-table summarize  \
 --i-table inptab.qza  \
 --o-visualization table.qzv \
 --m-sample-metadata-file ../../../mapping_LDA.txt
 
 #Collapse taxa levels at genus level

 qiime taxa collapse \
 --i-table inptab.qza \
 --i-taxonomy ../../features/taxonomy.qza \
 --p-level 6 \
 --o-collapsed-table inptab_genus.qza

module unload qiime2

date

Gneiss with OLS and LME

#!/bin/bash

# ----------------SLURM Parameters----------------

#SBATCH -J q2_gneiss
#SBATCH --mem=20g
#SBATCH --time=21:00:00
#SBATCH --ntasks=1
#SBATCH -n 8
#SBATCH -D /ufrc/strauss/andrea.nuzzo/projects/Workshop/merged
#SBATCH -o logs/5a_q2_gneiss_%j.out
#SBATCH -A strauss

date;hostname;pwd

################################################################################
#
# This script performs Gneiss log-transformed balances, then  OLS or LME modelling on your data
# 
################################################################################

# ----------------Housekeeping--------------------
cd models
rm -r gneiss
mkdir gneiss
cd gneiss

# ----------------Load Modules--------------------
module load qiime2

# ------------------Commands----------------------
 qiime feature-table filter-features \
   --i-table ../input/inptab.qza \
   --o-filtered-table filtered-table.qza \
   --p-min-samples 5

 qiime gneiss correlation-clustering \
   --i-table filtered-table.qza \
   --o-clustering hierarchy.qza

 qiime gneiss ilr-hierarchical \
   --i-table filtered-table.qza \
   --i-tree hierarchy.qza \
   --o-balances balances.qza

 qiime gneiss dendrogram-heatmap \
   --i-table filtered-table.qza \
   --i-tree hierarchy.qza \
   --m-metadata-file ../../../mappings_LDA.txt \
   --m-metadata-column Biostimulant \
   --p-color-map viridis \
   --o-visualization heatmap.qzv

 qiime gneiss ols-regression \
   --p-formula "TimePoint+Biostimulant+Biostimulant_TimePoint+Area+SPAD+Leaves+Stems+Root+Tot+LDA1_area_vs_treatment+LDA1_biomass_vs_treatment" \
   --i-table balances.qza \
   --i-tree hierarchy.qza \
   --m-metadata-file ../../../mappings_LDA.txt \
   --o-visualization regression_summary.qzv

 qiime gneiss lme-regression \
   --p-formula "Biostimulant+Area+SPAD+Leaves+Stems+Root+Tot+LDA1_area_vs_treatment+LDA1_biomass_vs_treatment" \
   --p-groups  TimePoint \
   --i-table balances.qza \
   --i-tree hierarchy.qza \
   --m-metadata-file ../../../mappings_LDA.txt \
   --o-visualization lme_all_vs_time.qzv
   
module unload qiime2

date

Exctract Balances

Option 1: one by one, qiime2 way

#!/bin/bash

# ----------------SLURM Parameters----------------

#SBATCH -J q2_balance_tax
#SBATCH --mem=10g
#SBATCH --time=1:00:00
#SBATCH --ntasks=1
#SBATCH -n 2
#SBATCH -D /ufrc/strauss/andrea.nuzzo/projects/Workshop/merged
#SBATCH -o logs/5b_balance_taxonomy_%j.out
#SBATCH -A strauss

date;hostname;pwd

################################################################################
#
# Finally, we want to extract the taxonomy for the balances which were 
# Significant in the previous steps
# 
################################################################################

# ----------------Housekeeping--------------------
cd models/gneiss
rm -r y*.qzv

# ----------------Load Modules--------------------
module load qiime2

# ------------------QUESTIONS---------------------
# Before going on, try to look at the summary of the model and answer: 
# How much is explained by the model? (R-squared)
# Does the model overfit? (pred_mse < model_mse)
# Do residuals show a trend? (if yes, you missed a variable)
# Are predicted in the residuals range? (Usually yes is expected for low R-squared)
# Do you have any significant balance in the OLS? (red in the heatmap)

# ------------------Commands----------------------
# If you have any significant balance and your model does not overfit, do the
# following for each balance to know who are the strains that significantly weigh
# for the chosen variable (i.e. Location)

qiime gneiss balance-taxonomy \
      --i-table composition.qza \
      --i-tree hierarchy.qza \
      --i-taxonomy ../../features/taxonomy.qza \
      --p-taxa-level 6 \
      --p-balance-name y0 \
      --m-metadata-file ../../../mappings_LDA.txt \
      --m-metadata-column TimePoint \
      --o-visualization y0_taxa_TimePoint.qzv 

qiime gneiss balance-taxonomy \
          --i-table composition.qza \
          --i-tree hierarchy.qza \
          --i-taxonomy ../../features/taxonomy.qza \
          --p-taxa-level 6 \
          --p-balance-name y0 \
          --m-metadata-file ../../../mappings_LDA.txt \
          --m-metadata-column Area \
          --o-visualization y0_taxa_Area.qzv
date

Option 2: Do a batch extraction with a python script and then analize the mastodontic csv file

I wrote a small python script that can be called through a bash job. As long as you kept the folder structure consistent with the Qiime2 tutorial it will work on your files as well. The script is hosted on GitHub so you can download and modify it if needed.

#!/bin/bash

# ----------------SLURM Parameters----------------

#SBATCH -J q2_gneiss
#SBATCH --mem=20g
#SBATCH --time=21:00:00
#SBATCH --ntasks=1
#SBATCH -n 8
#SBATCH -D /ufrc/strauss/andrea.nuzzo/projects/Workshop/
#SBATCH -o merged/logs/5c_gneiss_extractor_%j.out
#SBATCH -A strauss

date;hostname;pwd

################################################################################
#
# Extraction of all balances having adjusted p-values<0.01 from the gneiss viz
#
################################################################################

module load qiime2/2018.6

wget -O merged/scripts/Balance_extractor.py https://andreanuzzo.github.io/Strausslab/Balance_extractor.py
merged/scripts/Balance_extractor.py

module purge
date

To elaborate the mastodontic file you had to use either R again or Tableau or similar. Excel will fail poorly.

rel.abund <- transform_sample_counts(phylobj, function(x){x/sum(x)})
library(data.table)
## 
## Attaching package: 'data.table'
## The following objects are masked from 'package:dplyr':
## 
##     between, first, last
## The following object is masked from 'package:purrr':
## 
##     transpose
extracted_balances <- fread('merged/models/gneiss/ols_summary_dir/extracted_balances.csv')
Area.balances <-extracted_balances %>%
  dplyr::select(kingdom, phylum, class, order, family, genus, species) %>%
  mutate_at(vars(kingdom, phylum, class, order, family, genus, species), funs(gsub(".*__", "",.)))

subset_taxa(rel.abund, Genus %in% Area.balances$genus &
            !Kingdom =='Unassigned' & !Kingdom =='Archaea' & !Kingdom =='Plantae') %>%
  psmelt() %>%
  group_by(TimePoint, Biostimulant, Kingdom, Phylum, Class, Area) %>%
  summarize(mean=mean(Abundance, na.rm = TRUE)*100) %>%
  filter(mean>1e-2) %>%
  mutate(bin=ntile(Area,10)*round(min(Area, na.rm = TRUE),0)) %>%
  ggplot() +
      aes(x=Biostimulant, y = mean, colour=Class, size=bin) +
      geom_point(position = position_jitter()) +
      facet_grid(Kingdom~TimePoint) +
      labs(y = "Relative Abundance %", x = "Area (sqcm)") + 
      theme(axis.text.x=element_text(angle=90,hjust=1,vjust=0.5),
            legend.box = 'vertical') +
      guides(color=guide_legend(ncol=3)) +
      scale_colour_discrete(guide = guide_legend(title.position = "top", nrow = 1)) +
      scale_color_pander()+
      ggtitle('Classes differing for Area')
## Scale for 'colour' is already present. Adding another scale for
## 'colour', which will replace the existing scale.

Option 3: Jupyter notebook

The other option is to use the Qiime2 API and use a Jupyter notebook or an RStudio notebook to elaborate the balances. Requires a bit of python knowledge (you can do it!)

nano merged/scripts/jupyter_notebook.sh

And then

#!/bin/bash

#----------SLURM Parameters-----------------#

#SBATCH --job-name=jupyter_andy
#SBATCH -D /ufrc/strauss/andrea.nuzzo/Workshop/merged
#SBATCH --account=strauss
#SBATCH --output=logs/jupyter_%j.log
#SBATCH --ntasks=1
#SBATCH --cpus-per-task=4
#SBATCH --mem=10gb
#SBATCH --time=72:00:00

date;hostname;pwd

module load qiime2/2018.6
 
jupyter-notebook --no-browser --port=23456 --ip='*'
 
date
sbatch merged/scripts/jupyter_notebook.sh

Wait a couple of minutes. You will then have to open the log. You then need to open an ssh tunnel, using the address that the log gives you. It will appear as something like The Jupyter Notebook is running at:http://(c23a-s40.ufhpc or 127.0.0.1):23456/. Copy the part before ufhpc and complete this command in a new terminal window.

ssh -NL 8000:c21b-s18.ufhpc:23456 your_email@hpg.rc.ufl.edu

You will be prompted with the request for password. After that, open a new Firefox window and digit localhost:8000 to enter. If you didn’t set your password, use the token that it’s in the jupyter notebook log you open previously.

Random Forest Analyses

#!/bin/bash

# ----------------SLURM Parameters----------------

#SBATCH -J q2_RF
#SBATCH --mem=30g
#SBATCH --time=21:00:00
#SBATCH --ntasks=1
#SBATCH -n 8
#SBATCH -D /ufrc/strauss/andrea.nuzzo/projects/Workshop/merged
#SBATCH -o logs/6_q2_RF_%j.out
#SBATCH -A strauss

date;hostname;pwd

################################################################################
#
# This script performs Regression/classification with RF on your data
# 
################################################################################
# ----------------Housekeeping--------------------
cd models
rm -r randomforest
mkdir randomforest
cd randomforest

# ----------------Load Modules--------------------
module load qiime2/2018.6

# ------------------Commands----------------------

#Classification on categorigal variables

qiime sample-classifier classify-samples \
    --i-table ../input/inptab_genus.qza \
    --m-metadata-file ../../../mappings_LDA.txt \
    --m-metadata-column TimePoint \
    --p-optimize-feature-selection \
    --p-parameter-tuning \
    --p-estimator RandomForestClassifier \
    --p-n-estimators 500 \
    --p-cv 5 \
    --p-random-state 42 \
    --p-n-jobs -1 \
    --p-palette GreenBlue \
    --output-dir RFC

qiime sample-classifier regress-samples \
  --i-table ../input/inptab_genus.qza \
  --m-metadata-file ../../../mappings_LDA.txt \
  --m-metadata-column Area \
  --p-optimize-feature-selection \
  --p-parameter-tuning \
  --p-estimator RandomForestRegressor \
  --p-n-estimators 500 \
  --p-cv 5 \
  --p-random-state 42 \
  --p-n-jobs -1 \
  --output-dir RFR

module unload qiime2
date

Elaboration of random forest results

rel.abund <- transform_sample_counts(phylobj, function(x){x/sum(x)})
#These are the feature importances >1% in the RandomForest analysis made on qiime2
time.importance <-read.delim('merged/models/randomforest/Timepoint_feature_importance.tsv') %>%
  separate(feature, c("Kingdom", "Phylum", "Class", "Order", "Family","Genus"), sep=';') %>%
  mutate_at(vars(Kingdom, Phylum, Class, Order, Family, Genus), funs(gsub(".*__", "",.))) %>%
  filter(importance>1e-2)

subset_taxa(rel.abund, Genus %in% time.importance$Genus &
            !Kingdom =='Unassigned' & !Kingdom =='Archaea') %>%
  psmelt() %>%
  group_by(TimePoint, Biostimulant, Kingdom, Genus) %>%
  summarize(mean=mean(Abundance, na.rm = TRUE)*100, 
            stdev=sd(Abundance, na.rm = TRUE)) %>%
  ggplot() +
      aes(x=TimePoint, y = mean, fill=Genus) +
      geom_bar(stat = 'identity') +
      scale_fill_pander() + 
      facet_grid(Kingdom~Biostimulant) +
      theme(axis.text.x=element_text(angle=90,hjust=1,vjust=0.5)) +
      labs(y = "Relative Abundance %", x = "Filtered Samples") + 
  ggtitle('Species differing over time')

Area_RF <- read.delim('merged/models/randomforest/Area_feature_importance.tsv') %>%
  separate(feature, c("Kingdom", "Phylum", "Class", "Order", "Family","Genus"), sep=';') %>%
  mutate_at(vars(Kingdom, Phylum, Class, Order, Family, Genus), funs(gsub(".*__", "",.))) %>%
  filter(importance>1e-2)

subset_taxa(rel.abund, Genus %in% Area_RF$Genus & 
            !Kingdom =='Unassigned' & !Kingdom =='Archaea') %>%
  psmelt() %>%
  group_by(TimePoint, Area, Biostimulant, Kingdom, Genus) %>%
  summarize(mean=mean(Abundance, na.rm = TRUE), 
            stdev=sd(Abundance, na.rm = TRUE), 
            Area.mean = mean(Area, na.rm = TRUE),
            Area.stdev = sd(Area, na.rm = TRUE)) %>%
  write.csv('merged/models/randomforest/Area_RF.csv')